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Abstract The complex Cladonia mediterranea belongs to the section Impexae and is formed by C. azorica, C. maca- 
ronesica and C. mediterranea. These species are basically distributed in the Mediterranean and Macaronesian 
Regions. In the present work the limits between the species of this complex are re-examined. To this end, the mor- 
phological characters were studied along with the secondary metabolites and the DNA sequences from three loci 
(ITS rDNA, IGS rDNA and rpb2). The morphological data were studied by principal component analysis (PCA), while 
the DNA sequences were analyzed using several approaches available to delimit species: genealogical concordance 
phylogenetic species recognition, species tree (BEAST* and spedeSTEM) and cohesion species recognition. In 
addition, the genealogical sorting index was used in order to assess the monophyly of the species. The different 
procedures used in our study turned out to be highly congruent with respect to the limits they establish, but these 


taxonom 
y limits are not the ones separating the prior species. Either the morphological analysis or the different approaches 
to species delimitation indicate that C. mediterranea is a different species from C. macaronesica, while C. azorica 
and C. macaronesica, which are reduced to synonyms of C. portentosa, constitute a separate lineage. 
Article info Received: 20 February 2015; Accepted: 24 March 2015; Published: 30 April 2015. 
INTRODUCTION hypotheses for the evolutionary lineages to be considered as 


The development of molecular tools has brought about a more 
accurate species delimitation and a better understanding of the 
evolution of fungi. The definition of many species has changed. 
It is well-known that in many fungal groups a large number of 
morphological species hide cryptic species complexes (Bick- 
ford et al. 2007, Crespo & Lumbsch 2010). Despite the major 
methodological advances made in species delimitation, it still 
constitutes a challenge, since there does not exist a valid 
method that allows identification of independent evolutionary 
lineages in all the cases (Sites & Marshall 2003, Carstens et al. 
2013). One of the most used criterion for species delimitation 
in fungi has been the Genealogical Concordance Phylogenetic 
Species Recognition (Taylor et al. 2000), that uses several 
unlinked loci and reciprocal monophyly to identify the species. 
In many cases this criterion has been very useful to distinguish 
divergent lineages (Kroken & Taylor 2001, Dettman et al. 2003, 
Ott et al. 2004, Fournier et al. 2005, Doyle et al. 2013, Morgado 
et al. 2013). Nevertheless, due to species divergence being a 
temporal process, this criterion can fail in cases of delimitation 
of closely related species that have diverged recently (Knowles 
& Carstens 2007). Some other facts such as hybridization, re- 
combination and horizontal transfer can also cause the gene 
tree to be inconsistent with the species tree (Eckert & Carstens 
2008). There are of course other procedures used for species 
delimitation in fungi (Wirtz et al. 2008, Gazis et al. 2011). One 
of them is Templeton's (1989) cohesion species recognition, 
that does not require species monophyly (Weisrock & Larson 
2006, Wirtz et al. 2008, 2012). This method evaluates two 
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species. The first of them (H1) is that there is only one evolu- 
tionary lineage in the studied group; the second (H2) is that the 
evolutionary lineages are genetically or ecologically exchange- 
able. The rejection of both hypotheses along with the congru- 
ence of H2 with the lineages found in H1 permits to delimit the 
cohesion species (Templeton et al. 2000). Numerous methods 
based on coalescence have recently been combined with the 
species delimitation procedures (O'Meara et al. 2006, Pons et 
al. 2006, Liu et al. 2009, Ence & Carstens 2011, Yang & Ran- 
nala 2010). They have the advantage of taking into account 
the incomplete lineage sorting and not requiring the reciprocal 
monophyly (Carstens & Knowles 2007, Fujita et al. 2012). 
These methods have been already applied in several works 
on species delimitation in lichenized fungi (Leavitt et al. 2011, 
2012, 2013, Parnmen et al. 2012). 


An emergent approach that is gathering increasing approval is 
the one that uses diverse data and analysis types to trace the 
most significant evidence, which permits the establishing of 
boundaries among species (Padial & de la Riva 2010, Gazis et 
al. 2011, Gebiola et al. 2012). This is the procedure that some 
authors call ‘taxonomical integration’ (Wiens & Penkrot 2002, 
Dayrat 2005, Will et al. 2005, Padial et al. 2009). According 
to Carstens et al. (2013), several analytical methods must be 
used in order to delimit species within a group of organisms, 
since each of the extant methods takes prior assumptions 
that not always fit the available data or the speciation scenery 
under screening. 


Cladonia comprises most species within the family Cladonia- 
ceae. According to Stenroos et al. (2002), Cladonia is a mono- 
phyletic genus that encompasses all the species of the former 
genus Cladina (Ahti 1961, 1984, 2000), represented by about 
36 species from all over the world (Ahti 2000). This group, 
commonly known as reindeer lichens, is characterized by a 
crustose evanescent primary thallus, densely branched podetia, 
ecorticated, arachnoid surface, and by the absence of scyphi 
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and soredia (Ahti 1961, 1984). Furthermore, Stenroos et al. 
(2002) showed that Cladina is a monophyletic group (Group 
Cladinae) divided into two clades, one corresponding to the 
old section /mpexae and the other to the sections Cladina 
and Tenues. However, some studies indicated that the Group 
Cladinae is not monophyletic (DePriest et al. 1999, 2000, Guo 
& Kashiwadani 2004) but the old section /mpexae is monophy- 
letic. The section /mpexae is represented in Europe by C. azo- 
rica, C. macaronesica, C. mediterranea, C. portentosa and C. stel- 
laris. The problematic C. mediterranea complex includes three 
ofthese species, viz. C. azorica, C. macaronesica and C. medi- 
terranea. Cladonia azorica is reported to be widespread in Ma- 
deira, Azores, Ireland, England and Iceland (Ahti 1961, Ruoss 
1989, James 2009, Ahti & Stenroos 2013), while C. maca- 
ronesica is known from the Canary Islands, Madeira and 
Azores (Ahti 1961). Cladonia mediterranea has the broadest 
distribution, from Portugal to Turkey, southwestern Britain and 
the Canary Islands (des Abbayes & Duvigneaud 1947, Ruoss 
1989, James 2009, Hernández-Padrón & Pérez-Vargas 2009). 


It is still unclear whether C. azorica, C. macaronesica and C. medi- 
terranea represent independent species or not. Ahti (1977) syno- 
nymized C. azorica with C. macaronesica; but later on, he 
again recognized C. azorica (Ahti 1984), whereas Ruoss (1989) 
concluded that C. mediterranea and C. macaronesica were 
conspecific. However, in many floristic works C. macarone- 
sica and C. mediterranea are still treated as separate species 
(Hafellner 1995, Flores Rodrígues & Aptroot 2005, Carvalho 
et al. 2008, Hernández-Padrón & Pérez-Vargas 2009, Gabriel 
2012). Sicilia et al. (2009) refer to C. mediterranea group be- 
cause of the high morphological variation they found, while 
they pointed out the necessity of molecular studies for clarifying 
the taxonomy of this complex. To date, analyses using DNA 
sequence data have not been carried out. 


The aim of the present work is to study the species delimita- 
tion in the Cladonia mediterranea complex using different 
approaches and several data sources: DNA sequences from 
three loci (ITS rDNA, IGS rDNA and rpb2), morphological data 
and chemical data. 


MATERIAL AND METHODS 


Sampling 

The specimens were collected during 2011 from the Canary 
Islands, Madeira, Azores and the coast of Portugal. To com- 
plete the sampling, specimens deposited at MACB and H were 
included. In all, 40 specimens of C. azorica, C. macaronesica, 
C. mediterranea and C. portentosa were selected (Table 1). 
We included C. portentosa, which is a common species in the 
Iberian Peninsula and Macaronesia, because of its morphologi- 
cal resemblance to the other three species (Ahti 1961, Ruoss 
1989, Orange 1993, Burgaz & Ahti 2009, Sicilia et al. 2009, 
Ahti & Stenroos 2013). Two specimens of C. pycnoclada and 
two of C. confusa were also included, both South American 
members of the section /mpexae; they have been considered 
by some authors to be close to C. azorica (des Abbayes 1946, 
Tavares 1952). Cladonia deformis, C. boryi and C. cenotea were 
chosen as outgroup species, according to the results of Sten- 
roos et al. (2002). 


Morphology and secondary metabolites 


The morphological characters studied were selected on the 
basis of the original descriptions of the species (des Abbayes 
& Duvigneaud 1947, Ahti 1961, 1978) and other studies (Ruoss 
1989, Burgaz & Ahti 2009). Fourteen quantitative morphological 
characters were measured (length of podetia, width of pode- 
tia, number of branches, dichotomous branches percentage, 


trichotomous branches percentage, tetrachotomous branches 
percentage, branching angles, length of internodes, length of 
last branch, thickness of podetia, thickness of medulla, thick- 
ness of stereome, number of open axils, number of closed 
axils). For each specimen, the measures were taken in one to 
three podetia, according to the available material. All the mac- 
roscopic characters were measured by means of a digital slide 
gauge (0.01 mm precision) under a binocular stereomicroscope 
(Olympus SZX9). Transverse sections of the podetia were 
made free-hand, and the microscopical measures were taken 
at 400x magnification using an Olympus CX41 microscope, in 
distilled water. The matrix containing the fourteen characters 
previously mentioned was analyzed using principal component 
analysis (PCA). This analysis was performed with the Canoco 
4.5 program (ter Braak & Smilauer 2002). The variables were 
centered and standardized before the PCA. The values for the 
first two components were plotted (Fig. 1). In Fig. 1a the data 
were coded according to the morphological identification. Using 
the same matrix the discriminatory descriptors were inferred 
from the lenght of the vector and its correlation with the respec- 
tive axes, so Fig. 1b represents the correlation of the different 
morphological variables with the components. 


The secondary metabolites were studied in all the specimens 
using the solvents A (Toluene: dioxane: acetic acid) and B 
(Hexane: diethylether: formic acid) (White & James 1985). 


DNA extraction, PCR and sequencing 


Genomic DNA was extracted using DNeasy Plant Mini Kit 
(QIAGEN, Hilden, Germany), following the manufacturer's 
instructions. The DNA was eluted in the final step in 200 ml 
of elution buffer provided by the manufacturer. The following 
three nuclear loci were sequenced: ITS rDNA, rpb2 and IGS 
rDNA. The primers and PCR programs are described in Pino- 
Bodas et al. (2013). The amplifications were carried out with 
Ready-to-Go-PCR Beads (GE Healthcare Life Sciences, UK). 
PCR products were purified with ExoSAP-IT (USB Corpora- 
tion, OH, USA). The sequencing was performed at Macrogen 
Europe service (www.macrogen.com), with the same primers 
used for the PCR. 


Phylogenetic analyses 


The consensus sequences from forward and reverse templates 
were assembled and edited in Sequencher™ 4.1.4 program 
(Gene Codes Corporation, Inc, Ann Arbor, Michigan, USA). 
The sequences of each locus were manually aligned in BIO- 
EDIT 7.0 (Hall 1999). No ambiguous positions were found 
and all the positions of the alignments were included in the 
analyses. Each region was analyzed separately by Maximum 
Likelihood (ML) using RAxML 7.0.3 (Stamatakis et al. 2005), 
under the model GTRGAMMA. Fast bootstrap was run with 
500 pseudoreplicates. The congruence among the different 
topologies inferred from the loci was tested following Hillis & 
Bull (1993): each clade with more than 75 % bootstrap sup- 
port was scanned for conflict among loci. We considered the 
existence of a conflict whenever a clade was supported with a 
bootstrap (more than 75 %) in one locus, while it was not sup- 
ported in another locus, and the individual sequences of this 
clade were part of another clade with bootstrap support 75 96. 
In the combined datasets, only the specimens with sequences 
atleast for 2 genes were included. The combined dataset was 
analyzed by Maximum Parsimony (MP), ML and Bayesian 
analyses. MP analyses were performed in PAUP* v. 4.0.b.10 
(Swofford 2003) using the heuristic searches with 1 000 ran- 
dom taxon-addition replicates, with TBR branch swapping and 
MulTrees option in effect, equally weighted characters. Gaps 
were treated as missing data. For the confidence analysis, the 
bootstrap (Felsentein 1985) was applied, with 1 000 replicates 
and heuristic searchs. 
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a + C. azorica 


e C. mediterranea 


m C. portentosa 


^. C. macaronesica 


Fig. 1 Results of PCA analysis. a. PCA scatterplot with all the standardized variables and species studied; b. vector plot on PC1 and PC2. The length of 
vectors show the importance of each character. Bang: branching angles; Dich: dichotomous branches percentage; Labr: length of last branch; Lgth: length of 
podetia; Lint: length of internodes; Nbra: number of branches; Nucl: number of closed axils; Nuop: number of open axils; Tetr: tetrachotomous branches per- 
centage; Thme: thickness of medulla; Thpo: thickness of podetia; Thst: thickness of stereome; Tric: trichotomous branches percentage; With: width of podetia. 


The ML analysis was performed in the same conditions as the 
single gene datasets but considering 5 partitions: ITS rDNA, IGS 
rDNA and each codon position of rpb2. The best fit substitution 
model for each region was calculated using MrModeltest 2.3 
(Nylander 2004), with Akaike information criterion. The models 
selected and used in the Bayesian analysis were: GTR+G for 
IGS rDNA, SYM+G for ITS rDNA and SYM-I for rpb2. The 
Bayesian analysis was carried out using MrBayes 3.2 (Ronquist 
et al. 2012). The posterior probabilities were approximated by 
sampling trees using Markov Chain Monte Carlo (MCMC). The 
posterior probabilities of each branch were calculated by count- 
ing the frequency of trees visited during the MCMC analysis. 
Two simultaneous runs with 10 000 000 generations, each start- 
ing with a random tree and employing 4 simultaneous chains, 
were executed. Every 1 000th tree was saved into a file. The 
convergence was assessed checking that the average stand- 
ard deviation of split frequencies was < 0.01. In addition, the 
compare and slide commands in AWTY (Nylander et al. 2008) 
were used. Afterwards, the 50 % majority-rule consensus tree 
was calculated after removing the first 2 500 000 generations 
(i.e. the first 2 500 trees) using the ‘burn in’ command. 


Tests of monophyly and genealogical sorting index 


In order to assess the hypotheses that C. azorica, C. macaro- 
nesica and C. portentosa were monophyletic, constraint trees 
were constructed. These alternatives topologies were sup- 
plied to RAxML to search the ‘best’ ML tree. The GTRGAMMA 
model was used. Shimodaira-Hasegawa test (SH; Shimodaira 
& Hasegawa 1999) and expected likelihood weight test (ELW; 
Strimmer & Rambaut 2002) were performed using the program 
TREE-PUZZLE 5.2 (Schmidt et al. 2002) with 1 000 replicates 
resampled using the RELL method. 


The genealogical sorting index (GSI) was used to assess the 
level of genealogical exclusivity (Cummings et al. 2008) for 
C. azorica and C. macaronesica. The GSI value was not cal- 
culated for C. portentosa because we have more specimens of 
this species than of C. azorica and C. macaronesica, and the 
GSI has bias when unequal sampling size among the groups. 
The GSI was calculated for the ML tree of each locus and 
GSI, was calculated for the set of ML trees of the three loci. 


The significance was calculated using 10 000 permutations 
on the online platform at http://www.genealogicalsorting.org/. 
The GSI supports monophyly when this value is > 0.90 and the 
P-value < 0.001 according with Cummings et al. (2008) and 
Gazis et al. (2011). 


Species tree 


Two methods were used to calculate the species tree. First, 
the species tree under multispecies coalescent model was 
estimated using *BEAST implemented in BEAST (Drummond & 
Rambaut 2007), including only the specimens with sequence for 
the three loci. The specimens were assigned to species based 
on their morphology (C. azorica, C. macaronesica, C. mediter- 
ranea and C. portentosa). The model GTR+G was assigned to 
ITS and IGS partitions, and GTR+I to rpb2 partition, selecting 
birth-death speciation process with an uncorrelated relaxed 
lognormal clock (Drummond et al. 2006) and a constant 
size population. For the remaining priors the defect values 
were kept. The analysis was run for 50 000 000 generations, 
sampling every 1 000. The convergence was calculated with 
TRACER 1.5 (Rambaut & Drummond 2007). After discarding 
the first 10 000 000 generations the effective sample size for 
all the parameters of the evolutionary model reached values 
> 200. The tree was summarized with TREEANNOTATOR 1.7.5 
(Rambaut & Drummond 2013) using maximum clade credibility 
tree option for the target tree type. 


In the second method SpedeSTEM (Ence & Carstens 2011) 
was used to calculate the species tree. This method is based 
on coalescence that applies several loci gene trees to calculate 
the maximum likelihood species tree (Kubatko et al. 2009). This 
program allows not only to validate the species generated by 
other procedures, but also to delimit species with no a priori 
assignment of individuals. In this study only discovering analy- 
ses were made according to Satler et al. (2013). The ML gene 
trees were generated in PAUP (using the models estimated in 
MrModeltest), including C. mediterranea, C. macaronesica, C. por- 
tentosa, C. azorica and the outgroup. SpedeSTEM requires the 
specimens to be present in all the gene trees, and so only those 
can be studied for which it was possible to generate sequences 
for the three loci. SpedeSTEM needs a 9 value for scaling the 
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branch lengths in the species trees it produces. The 0 value 
for each locus was calculated in DnaSP v. 5 (Librado & Rozas 
2009), being 0 = 0.04437 for IGS rDNA, 6 = 0.04073 for ITS 
and 0 = 0.02804 for rpb2. The analysis was repeated for sev- 
eral 0 values: the average value of the three loci (8 = 0.03771), 
0 = 0.02, 0 = 0.03 and 0 = 0.04, to avoid the issues the program 
can have when calculating the likelinood for low 6 values (Giarla 
et al. 2014). 


Species delimitation by cohesion species recognition 


Haplotypes networks under statistical parsimony with a 
confidential interval of 95 % were generated with TCS 1.21 
(Clement et al. 2000) for each locus (ITS rDNA, IGS rDNA 
and rpb2). For the ITS rDNA analysis all the sequences of 
C. portentosa from GenBank (FR799166, FR799167, JQ695921, 
JQ695922, JQ695323) were included. Gaps were coded as 
missing data. The haplotypes were gathered manually in clades 
according to the rules of Templeton et al. (1987). This algorithm 
identifies clades united by mutational steps. The x-step clades 
are successively grouped in x+ 1-step clades and the final level 
of nested clades includes the complete network. The loops were 
resolved following the three criteria (frequency, topology and 
geographical) proposed by Pfenninger & Posada (2002). Table 
of contingence and Kruskal-Wallis analyses were done to test 
the null hypothesis (H2) of no significant association between 
haplotype variation and phenotypical variation. The quantitative 
variables with more contribution to separate the groups in the 
PCA analysis were analyzed by Kruskal-Wallis. The statistical 
analyses were performed in STATGRAPHICS 5.1 program. The 
clades 2-4 of ITS rDNA and 2-4 of rpb2 could not be included 
in the statistical analyses because they contained only one 
specimen. The pairwise fixation index F,, (Weir & Cockerham 
1984) was calculated with the program DnaSP. This value was 
used to assess whether gene flow exists among the cohesive 
species delimited. 


RESULTS 


Morphological analysis and secondary metabolites 


Fig. 1 shows the results of PCA. The first two principal compo- 
nents PC1 and PC2 summarize 51.44 % of the total variance 
(29.93 % and 21.51 96, respectively). The analysis distinguished 
two groups (Fig. 1a). The first group contains all the specimens 
of Cladonia mediterranea (on the upper left area of the scatter- 
plot) and the other group is formed by the rest of the species 
analyzed, C. azorica, C. macaronesica and C. portentosa (on 
the center of the scatterplot). The analysis shows a continuous 
morphological variation in the second group with a high degree 
of overlapping between the three species. The characters that 
most contribute to separate C. mediterranea from the other 
group were the dichotomous branching percentage and the 
number of closed axils (Fig. 1b). 


The secondary metabolites found in each specimen are listed 
in Table 1. All the specimens of C. mediterranea and C. maca- 
ronesica contained perlatolic and usnic acids. One specimen of 
C. azorica lacked usnic acid, containing only fumarprotocetraric 
and perlatolic acids. The other specimens contained fumar- 
protocetraric, perlatolic and usnic acids. Three specimens of 
C. portentosa lacked usnic acid (C. portentosa subsp. porten- 
tosa f. subimpexae). The other specimens contained perlatolic 
and usnic acids. 


Datasets and phylogenetic analyses 


In this study 113 new sequences (39 from ITS rDNA, 40 from 
IGS rDNA and 34 from rpb2) have been generated, the Gen- 
Bank accession numbers are listed in Table 1. The concate- 


nated dataset contained 1 781 characters, 136 of which were 
parsimony-informative. The locus which contained more 
parsimony-informative positions was ITS rDNA with 50 posi- 
tions, followed by IGS rDNA with 47 positions and rpb2 with 
39 positions. 


The separated analyses of ITS rDNA and IGS rDNA yielded 
trees with the same topology (not shown). Three well supported 
clades appeared. One clade contained all the samples of 
C. mediterranea, the second clade contained the samples of 
C. confusa and C. pycnoclada and the third included all the 
samples of C. azorica, C. macaronesica and C. portentosa. The 
ML analysis of rpb2 yielded a tree with 4 clades, one with all the 
samples of C. mediterranea, two clades containing samples of 
C. azorica, C. macaronesica and C. portentosa, and the fourth with 
the samples of C. confusa and C. pycnoclada. Only the C. con- 
fusa and C. pycnoclada clade was supported. No conflict 
among the loci was found and the datasets were combined. 
The MP analysis based on the concatenated dataset yielded 
trees of 357 steps long, Cl = 0.8711 and RC = 0.9170. The ML 
analysis yielded a tree with a likelihood value of -LnL = 4431.65, 
while the mean likelihood of the Bayesian tree sampling was 
-LnL = 4620.16. The tree topology was the same in all the 
analyses, so only the 50 % consensus majority tree from 
Bayesian analysis is shown (Fig. 2). Three clades appeared in 
all the analyses. One of them contained all the specimens of 
C. mediterranea, the second the specimens of C. pycnoclada 
and C. confusa, and the third the specimens of C. azorica, 
C. macaronesica and C. portentosa. The C. mediterranea 
clade was basal to the group. This clade was well supported 
in MP and ML analyses but not in the Bayesian analysis (0.74 
posterior probability). The clade formed by C. pycnoclada and 
C. confusa was supported in all the analyses, and the third clade 
had high support in all the analyses. In the ML and Bayesian 
analyses two unsupported subclades can be distinguished in 
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Fig.2 Phylogeny of Cladonia mediterranea complex. The 50 % consensus 
majority tree of the Bayesian analysis based on concatenated dataset. The 
values on the branches are the posterior probability from Bayesian analysis 
(2 0.95), bootstrap values from ML analysis (2 75 96) and bootstrap values 
from MP analysis (2 75 96). 
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Fig. 3 a. Species tree inferred in *BEAST based on the concatenated 
dataset. The values on the branches correspond to the posterior probability; 
b. species tree estimated in spedeSTEM based on discovery approach. 


C.azorica —— 
C. macaronesica 
1 C. portentosa 


C. macaronesica 


0.75 


the third clade. One of them was constituted by some samples 
of C. macaronesica and C. azorica and the other was formed 
by all the specimens of C. portentosa, some specimens of 
C. macaronesica and some specimens of C. azorica. 


Species trees 


The results of the *BEAST analysis are shown in Fig. 3a. Cla- 
donia mediterranea was well supported, C. azorica and C. por- 
tentosa form a clade but the node was not significantly sup- 
ported with posterior probability. The clade clustering C. azorica, 
C. macaronesica and C. portentosa was significantly supported. 
The results from SpedeSTEM analyses were similar for different 
0 values. In all the cases the model that obtained most support 
was the one that considers 2 lineages in the ingroup (Fig. 3b). 
One of them was composed only by C. mediterranea and the 
other by C. macaronesica, C. portentosa and C. azorica. The 
probability of this model was wi = 0.99 for 0 = 0.03771, wi = 1.0 
for 0 = 0.02, wi = 0.99 for 0 = 0.03 and wi = 0.99 for 0 = 0.04. The 
probability for alternative models was wi 7 0.0. 


Hypotheses and GSI 


The SH and ELW significantly rejected the three hypotheses: 
a. the monophyly of C. azorica (SH, P-values = 0.0090 and 
ELW P-value = 0.0009); 
b. the monophyly of C. macaronesica (SH, P-values = 0.0290 
and ELW P-value = 0.0423); 
c. the monophyly of C. portentosa (SH, P-value = 0.0270 and 
ELW P-value = 0.0142). 


The GSI test results are shown in the Table 2. The GSI values 
for C. azorica were similar among the different loci and the 
P-values rejected the monophyly in all the loci. The GSI values 
for C. macaronesica were low in ITS rDNA and rpb2, and not 
significant. However, the GSI value of IGS rDNA was 0.5806 
and significant. The GSI, rejected the exclusive ancestry for 
both species. 


Networks and nested clade analyses 


A total of fifteen haplotypes of ITS rDNA were identified, con- 
nected in a single network (Fig 4a). Two haplotypes were unique 
for C. mediterranea, two were unique for C. macaronesica, two 
were unique for C. azorica and five were unique for C. por- 
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tentosa. The other four haplotypes were shared by samples of 
different species (C. macaronesica, C. portentosa and C. azo- 
rica). The IGS rDNA network analysis yielded a total of six haplo- 
types connected into a single network (Fig. 4b). All the samples 
of C. mediterranea were represented in one unique haplotype, 
one haplotype was unique for C. azorica, and one was unique 
for C. portentosa. The other three haplotypes were shared by 
samples of C. macaronesica and C. azorica; C. macaronesica 
and C. portentosa; or C. macaronesica, C. portentosa and 
C. azorica. The rpb2 network analysis yielded nine haplotypes 
connected into a single network (Fig. 4c), four of them were 
unique for C. mediterranea, one for C. macaronesica and one 
was unique for C. azorica. The other three haplotypes were 
shared by samples of different species. 


The nested clade analyses generated five 1-step clades, four 
2-step clades and two 3-step clades for ITS rDNA; for IGS rDNA, 
three 1-step clades and two 2-step clades were generated; and 
for rpb2 six 1-step clades, four 2-step clades and two 3-step 
clades were generated. All the specimens of 2-2 clade from IGS 
rDNA and 3-1 clade from ITS were identified as C. mediterra- 
nea, while the 3-2 clade of rpb2 contained all the specimens of 
C. mediterranea and one of C. macaronesica. The specimens 
of C. macaronesica, C. portentosa and C. azorica appeared 
together in the 3-2 clade of ITS rDNA, 2-1 clade of IGS rDNA 
and 3-1 clade of rpb2. The specimens grouped together in the 
2-1 clade of IGS rDNA, 3-2 clade of ITS rDNA and 3-1 clade of 
rpb2 are from Macaronesia, North America and Europe while 
the specimens of the clades 2-2 of IGS rDNA and 3-2 of rpb2 
are from the Canary Islands and the Iberian Peninsula. 


Table 3 and 4 summarize the results of the contingency table 
analyses. In the analyses of the 3 networks (ITS rDNA, IGS 
rDNA and rpb2) significant differences in characters were ob- 
served (Table 3). These characters include medulla (loose/com- 
pact), branching pattern (isotomic/anisotomic/subisotomic) 
and algal layer (continuous/discontinuous). These differences 
were observed between the 3-step clades of ITS rDNA and 
rpb2 and 2-step clades in IGS rDNA. No significant differences 
in the presence/absence of fumarprotocetraric acid among 
these clades were found. The results of the contingency table 
analyses at 2-step level are presented in Table 4. In ITS rDNA, 
most of the significant differences were detected among the 
2-1 clade and the other clades; in rpb2 significant differences 


Table 2 Genealogical sorting index and probability values under the null hypothesis that the samples labeled as putative species are monophyletic. 


. ITS rDNA IGS rDNA rpb2 GSI, 
Species 
GSI P-value GSI P-value GSI P-value GSIt P-value 
C. azorica 0.1429 0.1798 0.2114 0.0659 0.2614 0.0463 0.1429 0.1782 
C. macaronesica 0.1556 0.1376 0.5806 1e-04* 0.1795 0.3308 0.1556 0.1343 


* denotes significant result. 
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Fig.4 Haplotype networks at 95 % of probability of C. mediterranea complex based on a. ITS rDNA; b. IGS rDNA; c. rpb2. The circles represent the haplotypes 
and the size is proportional to haplotype frequency. The small circles represent missing haplotypes. The discontinuous lines outline the 1-step nested clades, 
the grey lines outline the 2-step nested clades and the black lines outline the 3-step nested clades. 


were detected among the 2-3 clade and the other clades; and 
in IGS rDNA significant differences were detected among the 
1-3 clade and the other clades. 


Table 5 shows the Kruskal-Wallis results. Significant differences 
were obtained for all of the characters among the 3-step clades 
in ITS rDNA and rpb2 and the 2-step clades in IGS rDNA. How- 
ever, there were not significant differences among all the 2-step 
clades (see Tukey test, Table 6). No significant differences were 
found among the 2-step clades in rpb2 for dichotomous branch- 
ing rate and trichotomous branching rate. 

The F „values between the 3-step clades of ITS rDNA and rpb2 
and the 2-step clades of IGS rDNA are shown in Table 7. In all 
the comparations the values were high. The lowest value was 
between the clades appearing in rpb2. 


DISCUSSION 


This work addresses the species delimitation in the C. medi- 
terranea complex using two data sources: phenotypical data 
(morphological and chemical) and DNA sequences from three 
nuclear genes. The DNA data were analyzed by different meth- 
ods often used for species delimitation (gene trees, species 
trees, haplotype networks) and they were highly congruent. 
When the analyses performed with different type of data show 
congruent results (as in this case), the concordant inferred 
boundaries likely correspond to existing biological entities. Ac- 
cording to our results, the most probable scenario is the one 
that comprises two species. 


The analyses of the morphological data reveal that C. mediter- 
ranea is different from the remaining species. The genealogical 
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Table 3 Contingency test results for association of haplotype clades and 
phenotypical characters at highest step clade level for each locus. 


Table 5 Kruskal-Wallis results for association of haplotype clades and pheno- 
typical characters. 


Comparation level Character P-value Comparation level Character Estatistic P-value 
ITS rDNA ITS rDNA 
3-step Presence / lack of fumarprotocetraric acid 0.193 3-step Dichotomous branching rate (%) 5.53597 0.0186* 
3-step Isotomic / anisotomic or subisotomic pattern 0.001* 3-step Trichotomous branching rate (76) 3.90411 — 0.0481* 
3-step Compact / loose medulla 0.000* 3-step Closed axils rate (96) 12.4813  0.0004* 
3-step Continuous / discontinuous algal layer 0.000* 2-step Dichotomous branching rate (76) 8.26747  0.0407* 
IGS rDNA 2-step Trichotomous branching rate (%) 9.39649 0:0244" 
2-step Presence / lack of fumarprotocetraric acid 0.211 isis EIOS GONE tates) 12:9191 USD 
2-step Isotomic / anisotomic or subisotomic pattern 0.007* IGS rDNA 
2-step Compact / loose medulla 0.000* 2-step Dichotomous branching rate (96) 10.6712 . 0.0010* 
2-step Continuous / discontinuous algal layer 0.000* 2-step Trichotomous branching rate (76) 7.5518 0.0059* 
2-step Closed axils rate (96) 11.8636 . 0.0057* 
wpe j , 1-step Dichotomous branching rate (96) 10.9286 0.0042* 
dns Gieres umarpiotoceiaric:acid pee 1-step Trichotomous branching rate (%) ^ 7.64197 — 0.0219* 
3-step Isotomic / anisotomic or subisotomic pattern 0.004 1-step Closed axils rate (26) 43.5391  0.0014* 
3-step Compact / loose medulla 0.000* 
3-step Continuous / discontinuous algal layer 0.000* rpb2 
* denotes significant result. 3-step Dichotomous branching rate (96) 5.16176 . 0.02308* 
3-step Trichotomous branching rate (96) 4.64772  0.03109* 
3-step Closed axils rate (96) 11.268 0.00078* 
Table 4 Contingency test results for association of haplotype clades and 2-step Dichotomous branching rate (96) 5.63538 0.1307 
phenotypical characters. 2-step Trichotomous branching rate (96) 6.65459 0.0837 
2-step Closed axils rate (96) 14.1365 . 0.00272* 


Comparation level Character P-value 
ITS rDNA 

2-1 to 2-2 Presence / lack of fumarprotocetraric acid 0.429 
2-1 to 2-3 Presence / lack of fumarprotocetraric acid 0.171 
2-2 to 2-3 Presence / lack of fumarprotocetraric acid 0.412 
2-1 to 2-2 Compact / loose medulla 0.001* 
2-1 to 2-3 Compact / loose medulla 0.000* 
2-2 to 2-3 Compact / loose medulla 0.704 
2-1 to 2-2 Continuous / discontinuous algal layer 0.007* 
2-1 to 2-3 Continuous / discontinuous algal layer 0.000* 
2-2 to 2-3 Continuous / discontinuous algal layer 0.296 
2-1 to 2-2 Isotomic / anisotomic or subisotomic pattern 0.001* 
2-1 to 2-3 Isotomic / anisotomic or subisotomic pattern 0.002* 
2-2 to 2-3 Isotomic / anisotomic or subisotomic pattern 0.406 
IGS rDNA 

1-1 to 1-2 Presence / lack of fumarprotocetraric acid 0.179 
1-1 to 1-3 Presence / lack of fumarprotocetraric acid 0.091 
1-2 to 1-3 Presence / lack of fumarprotocetraric acid 0.348 
1-1 to 1-2 Compact / loose medulla 0.662 
1-1 to 1-3 Compact / loose medulla 0.002* 
1-2 to 1-3 Compact / loose medulla 0.487 
1-1 to 1-2 Continuous / discontinuous algal layer 0.450 
1-1 to 1-3 Continuous / discontinuous algal layer 0.018* 
1-2 to 1-3 Continuous / discontinuous algal layer 0.000* 
1-1 to 1-2 Isotomic / anisotomic or subisotomic pattern 0.134 
1-1 to 1-3 Isotomic / anisotomic or subisotomic pattern 0.160 
1-2 to 1-3 Isotomic / anisotomic or subisotomic pattern 0.001* 
rpb2 

2-1 to 2-2 Presence / lack of fumarprotocetraric acid 0.166 
2-1 to 2-3 Presence / lack of fumarprotocetraric acid 0.054 
2-2 to 2-3 Presence / lack of fumarprotocetraric acid 1.000 
2-1 to 2-2 Compact / loose medulla 0.045* 
2-1 to 2-3 Compact / loose medulla 0.000* 
2-2 to 2-3 Compact / loose medulla 0.022* 
2-1 to 2-2 Continuous / discontinuous algal layer 0.191 
2-1 to 2-3 Continuous / discontinuous algal layer 0.000* 
2-2 to 2-3 Continuous / discontinuous algal layer 0.006* 
2-1 to 2-2 Isotomic / anisotomic or subisotomic pattern 0.208 
2-1 to 2-3 Isotomic / anisotomic or subisotomic pattern 0.015* 
2-2 to 2-3 Isotomic / anisotomic or subisotomic pattern 0.003* 


* denotes significant result. 


phylogenetic species recognition (GPSR) was congruent with 
the results of the analysis of the morphological data (Fig. 2). 
Cladonia mediterranea formed a monophyletic clade well sup- 
ported in MP and ML analyses, but not in the Bayesian analysis. 
The hypotheses tests (SH and EWL) significantly rejected the 
alternative topologies, in which C. azorica, C. macaronesica 
and C. portentosa were divided into independent monophyletic 
groups. Since the incomplete lineage sorting could be responsi- 
ble for the lack of monophyly of C. azorica, C. macaronesica and 


* denotes significant result. 


Table 6 Tukey’s multiple comparison test for significant results of the 
Kruskal-Wallis analyses. 


Dichotomic Trichotomic Closed axil 
ITS rDNA 
2-1 to 2-2 i ns 7 
2-1 to 2-3 ns ns i 
2-1 to 2-4 ns ns ns 
2-2 to 2-3 ns ns ns 
2-2 to 2-4 ns ns ns 
2-3 to 2-4 ns ns ns 
IGS rDNA 
1-1 to 1-2 ns ns ns 
1-1 to 1-3 $ Y * 
1-2 to 1-3 ns ns ns 
rpb2 
2-1 to 2-2 - - ns 
2-1 to 2-3 - - ns 
2-1 to 2-4 - - * 
2-2 to 2-3 - - * 
2-2 to 2-4 - - ns 
2-3 to 2-4 - - ns 


ns = not significant; * = significant with 95 % of probability. 


Table7 Pairwise F,, for each clade defined in the networks. 


Locus Comparisons Es 

ITS rDNA 3-1 to 3-2 step clade 0.87659 
IGS rDNA 2-1 to 2-2 step clade 0.93114 
rpb2 3-1 to 3-2 step clade 0.69796 


C. portentosa, we applied the GSI test to evaluate the degree 
of genealogic divergence. The monophyly of C. azorica and 
C. macaronesica was not supported by this test. The species 
trees generated by means of *BEAST and SpedeSTEM gave 
rise to two well-supported species (Fig. 3). These analyses are 
congruent with the gene trees and the morphological analysis, 
leading to consider C. azorica, C. macaronesica and C. porten- 
tosa as a unique species, and C. mediterranea as a different 
one. The cohesion species recognition requires, in addition 
to rejecting the two null hypothesis, that the groups delimited 
during the evaluation of H1 be congruent with H2 hypothesis 
(Templeton et al. 2000). This congruence happens at 3-step 
clade level in ITS rDNA and rpb2 and at 2-step clade level in 
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IGS, since at an inferior level (2-step clade level in ITS and rpb2 
and 1-step clade level in IGS) significant results were obtained, 
but not among all the clades. The morphological differences 
occur between the clades that contain samples of C. mediter- 
ranea and the remaining clades, while there are no significant 
differences between the clades that contain the samples of 
C. azorica, C. macaronesica and C. portentosa (2-2, 2-3 and 
2-4 in ITS rDNA, 2-1, 2-2, 2-4 in rpb2 and 1-1 and 1-2 in IGS 
rDNA). Strong evidence for the fact that C. mediterranea is a 
different species from C. macaronesica, is that all the samples 
of C. mediterranea are confined to a unique clade in all the 
haplotype networks. In addition, C. mediterranea shows high 
levels of genetic differentiations, according with the F... values. 


The analyses of the morphological data and also numerous 
analyses based on the DNA sequences are consistent, indicat- 
ing that C. mediterranea is an independent evolutionary line- 
age and C. azorica, C. macaronesica and C. portentosa are 
conspecific. Thus our results reject the taxonomical proposal 
that C. mediterranea and C. macaronesica are conspecific 
(Ruoss 1989). This author studied the branching pattern and the 
characteristics of the algal layer and concluded that C. maca- 
ronesica and C. mediterranea were the same species. The di- 
agnostic characters used to distinguish these species were the 
following: length of the internodes are longer in C. mediterranea 
than in C. macaronesica; the algal layer is continuous in C. me- 
diterranea and discontinuous in C. macaronesica; a compact 
medulla present in C. mediterranea and a lax medulla in C. mac- 
aronesica; the axils are frequently closed in C. mediterranea 
and generally open in C. macaronesica (Ahti 1961). The PCA 
analyses carried out in this work show that the most relevant 
variables to distinguish C. mediterranea from the remaining 
species are the percentage of dichotomous branching and the 
number of closed axils. According to Ruoss (1989) C. me- 
diterranea had more closed axils than C. macaronesica. How- 
ever, we think that the internodal length does not contribute to 
the separation of C. mediterranea from the remaining samples, 
since C. portentosa has internodes of similar length or even 
with greater variation (Burgaz & Ahti 2009). Burgaz & Martínez 
(2008) found that the podetial wall is thicker in C. mediterranea 
than in C. portentosa; however, in our analysis this character 
had a scant contribution to distinguish C. mediterranea from the 
other species. The morphological characters that distinguish 
C. mediterranea are: the presence of a continuous algal layer, 
the presence of a compact medulla, the prevalence of isotomy, 
with dichotomous branching and closed axils. These characters 
are the originally used ones to describe the species (des Ab- 
bayes & Duvigneaud 1947). 


The boundaries among C. azorica, C. macaronesica and C. por- 
tentosa were not supported by any of the analyses carried out 
in this work. Cladonia azorica was distinguised from C. macaro- 
nesica mainly by the presence of fumarprotocetraric acid and by 
having a greater number of trichotomous branching, although 
the dichotomous pattern is also prevailing in this species (Ahti 
1961). But our analyses did not show a correlation between the 
presence of fumarprotocetraric acid and a greater number of 
trichotomous branching. In previous works based on the study 
of the morphological variation, the species status of C. azorica 
and C. macaronesica had already been questioned (Ahti 1977). 
Although no previous study has suggested that C. portentosa 
is conspecific with the latter, Orange (1993) pointed out that in 
Britain it was impossible to distinguish C. azorica from C. por- 
tentosa only by means of morphological characters. The mor- 
phological similarities of C. macaronesica and C. azorica with 
C. portentosa are clear, even with C. mediterranea (Fig. 5), and 
the identification keys and floristic works usually point out the 
possible confusion of C. portentosa with these other species 


(James 2009, Sicilia et al. 2009). But C. portentosa is generally 
distinguished by the prevailing trichotomous branching and an 
anisotomous pattern, where a main axis is clear. Nonetheless, 
C. portentosa is a very variable species, either morphologically 
or chemically (des Abbayes 1939, Ahti 1961, 1978, Burgaz & 
Martínez 2008). Within this taxon several forms and subspe- 
cies have been described. Cladonia portentosa subsp. pacifica, 
growing in western North America, is more slender and deflexed 
than C. portentosa subsp. portentosa and shows a greater 
number of dichotomous branches (Ahti 1978, Brodo & Ahti 
1996). Cladonia portentosa subsp. pacifica f. decolorans is a 
chemotype that lacks usnic acid, turning to a greyish shade 
(Brodo & Ahti 1996). Cladonia portentosa subsp. portentosa 
f. subimpexa also lacks usnic acid (Ahti 1978). The ITS rDNA 
sequences of C. portentosa subsp. pacifica and C. porten- 
tosa subsp. portentosa were recently compared and it was 
found that there was no genetical difference between them 
(Smith et al. 2012). Our analyses confirm these results. The 
two specimens of C. portentosa subsp. pacifica here included 
share a haplotype with some of C. portentosa subsp. portentosa 
samples in each of the 3 loci. 


In other species of the Group Cladinae (Stenroos et al. 2002) 
similar results have been found. This is the case of C. arbuscula, 
for which several subspecies were defined on the basis of the 
morphological and chemical variation. However, much of this 
variation is not correlated with the genetic variation (Piercey- 
Normore et al. 2010). The authors attribute the high variation 
within this species to environmental agents such as lighting, 
humidity, nutrients and thallus age. The warm temperatures 
throughout the year in Macaronesia, which causes a continu- 
ous development of the podetia, could be the environmental 
agent that determines C. portentosa to develop a prevailing 
dichotomous branching, instead of trichotomous. Ahti (1961) 
had already pointed out that in southern Europe (Portugal) 
C. portentosa tended to produce dichotomous branching, being 
easily mistaken for C. mediterranea, with which it often coexists 
(Burgaz & Ahti 2009). 


Our results indicate that C. confusa and C. pycnoclada are re- 
lated, while in the phylogeny submitted by Stenroos et al. (2002) 
C. confusa appeared closely related to C. terra-novae and 
C. portentosa. However, in our analyses C. confusa is not mono- 
phyletic, which could reveal a lack of genetic homogeneity of this 
species. Further studies, based on a wide range of sampling, 
should be made to confirm this observation. 


TAXONOMY 


In this section we present formally the taxonomical changes. 


Cladonia portentosa (Dufour) Coem., Bull. Acad, Roy. Sci. 
Belgique, ser. 2, 19: 49. 1865. 


Basionym. Cenomyce portentosa Dufour, Ann. Gen. Sci. Phys. 8: 69. 1821. 

Type. France, Landes, Saint-Sever-sur-Adour, 1818, J.-M. Dufour (PC- 
Herb. Desmaziéres lectotype, Ahti, Ann. Bot. Fenn. 15: 8. 1978). 

= Cladonia azorica Ahti, Ann. Bot. Soc. Zool.-Bot. Fenn. Vanamo 32, 1: 
36. 1961. 

Type. PoRrUcaL, Azores, Pico, 25 Sept. 1961, A.G. da Cunha & L. Sobrinho 
(LISU holotype). 

= Cladonia macaronesica Ahti, Ann. Bot. Soc. Zool.-Bot. Fenn. Vanamo 
32, 1: 37. 1961. 

Type. PoRrucA., Madeira, Entre as Queimadas e o Caldeirão Verde, 1951, 
C.N. Tavares 4583 (LISU holotype). 
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Fig. 5 Photographs of the four species studied, showing the general configurations of podetia a. Cladonia azorica (Haikonen 26865, H); b. C. macaronesica 
(Pérez-Vargas s.n., TFC 10602); c. C. mediterranea (Burgaz s.n., MACB 61559); d. C. portentosa (Stenroos 6074, H). 
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